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ABSTRACT 

Measurements of the intergalactic medium (IGM) temperature provide a potentially 
powerful constraint on the reionisation history due to the thermal imprint left by the 
photo-ionisation of neutral hydrogen. However, until recently IGM temperature mea- 
surements were limited to redshifts 2 < z < 4.8, restricting the ability of these data to 
probe the reionisation history at z > 6. In this work, we use recent measurements of the 
IGM temperature in the near-zones of seven quasars at z ~ 5.8 — 6.4, combined with 
a semi-numerical model for inhomogeneous reionisation, to establish new constraints 
on the redshift at which hydrogen reionisation completed. We calibrate the model to 
reproduce observational constraints on the electron scattering optical depth and the 
H i photo-ionisation rate, and compute the resulting spatially inhomogeneous temper- 
ature distribution at z ~ 6 for a variety of reionisation scenarios. Under standard 
assumptions for the ionising spectra of population-II sources, the near-zone tempera- 
ture measurements constrain the redshift by which hydrogen reionisation was complete 
to be z r > 7.9 (6.5) at 68 (95) per cent confidence. We conclude that future temper- 
ature measurements around other high redshift quasars will significantly increase the 
power of this technique, enabling these results to be tightened and generalised. 
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1 INTRODUCTION 

The reionisation of cosmic hydrogen represents a signifi- 
cant moment in the history of the Universe. The appear- 
ance of the first stars and quasars heated and ionised the 
IGM, catalysing the transition in which intergalactic hydro- 
gen changed from being predominantly neutral to its present 
day, highly ionised state. Understanding exactly how and 
when this phase transition occurred will offer insight into 



the nature of the first sources of light (e.g. Barkana & Loeb 
|2001[ ). However, existing constraints on the timescale and 
extent of reionisation are currently uncertain (Choudhury 



fc Ferrara|2006| |Pritchard et al.|2009| |Mitra et al.|2011 1 

There are two primary observations which are used to 
constrain the epoch of hydrogen reionisation. The first of 
these are the absorption spectra of quasars at z ~ 6 ( White 
et al.||2003| |Fan et al.||2006| |Becker et aT||2007[ ). These en- 
able detailed study of the IGM ionisation state at high red- 
shift, and indicate that the Universe was reionised at z > 6 
( Gnedin fc Fan|2006| but see also |Mesinger|2010| . The sec- 
ond observation is provided by measurements of the Thom- 
son optical depth, tomb, to the surface of last scattering, 



which indicate that hydrogen reionisation could not have 
ended any earlier than z r = 10.6±1.2 ( Komatsu et al.|2011 1. 
In practice, this constraint is consistent with a wide variety 
of extended reionisation scenarios. 

Analytical and semi-numerical models have been con- 
fronted with a combination of these and other observations 
to place more quantitative limits on the reionisation red- 



shift. In particular, studies by Choudhury & Ferrara ( 2006 1, 
|Pritchard et al.| ( |2009[ ) and |Mitra et al.f ( |2011| ) employ the 
above observables, alongside the observed distribution of Ly- 
man limit systems and constraints on the star formation 
rate, to model reionisation. These models can be used to in- 
fer constraints on the extent and completion time of reion- 
isation. However, the lack of direct observational evidence 
at z > 6 means that parameters in these models, such as 
the ionising background intensity and the thermal state of 
the IGM, are relatively uncertain at high redshift. This im- 
plies that, from current observations alone, constraints on 
the redshift by which reionisation is complete remain loose, 



and cannot be improved beyond 6 < z r < 10.5 (Mitra et al. 
[20TH >. 

A promising additional approach comes from measure- 
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ments of the IGM temperature ( Zaldarriaga et al. 2001 Trac 



et al. 2008 Furlanetto & Oh 2009 Bolton et al. 20101 



The relatively long time scale of adiabatic cooling implies 
that the low redshift IGM still carries the thermal imprint 
of reionisation heating from much earlier times. Several au- 
thors have previously used this fact, combined with temper- 
ature measurements at z < 4, to constrain the end of the 
hydrogen reionisation epoch to z r < 10 ( Theuns et al.|2002 
Hui fc Haiman|2 003). Unfortunately, because the imprint of 
reionisation heating at z r > 10 reaches a thermal asymp- 
tote by z ~ 4 these studies are limited by the relatively low 
redshift at which the temperature measurements were made. 

Measurements of IGM temperature at higher redshift 
therefore offer a new opportunity to tighten these con- 
straints. A first step towards this was made by |Bolton et al. | 
[] ( |2010[ ), who obtained a measurement of the IGM tempera- 
ture at mean density from Lya absorption lines in the near- 
zone of the z = 6.02 quasar SDSS J0818+1722. Although the 
statistical uncertainty on their single line-of-sight measure- 



ment was large, Bolton et al. ( 2010 1 were nevertheless able 



to use a simple cooling argument to infer a reionisation red- 
shift of z T < 9 (to 68 percent confidence) within the quasar 
near-zone, assuming hydrogen was reionised by population- 
II stars. Recently, Bolton et al. (20111 have significantly 



extended these observational constraints by measuring the 
IGM temperature in the near-zones of six further quasars 
at 5.8 < 2 < 6.4. In this work, we therefore seek to extend 
the modelling of |Bolton et al. ] ( |2010[ ) and include all seven 
near-zone temperature measurements to obtain a limit on 
the reionisation redshift for the general IGM. Specifically, 
by using a semi-numerical model to follow inhomogeneous 
reionisation and photo-heating of intergalactic hydrogen by 
stellar sources, as well any subsequent photo-heating by the 
quasar itself, we may calculate near-zone temperatures as a 
function of spatial position and reionisation redshift. A com- 
parison of our model with the observational data then allows 
us to constrain the redshift at which reionisation completed. 

The paper is organised as follows. We begin in Section[2] 
by introducing our semi-numerical model of reionisation, 
and proceed in Section [3] to give a more detailed description 
of the modelling of ionising sources. Predicted temperatures 
from the model for different reionisation scenarios are then 
described in Section [4] and we compare these to the [Bolton] 
et al. | ( [2011 1 quasar temperature measurements. Our con- 
straints on reionisation history are then discussed in Sec- 
tion [5j where we also assess the prospects for improving on 
these constraints in the future using additional quasar sight- 
lines. Conclusions and further directions are then presented 
in Section [6] Throughout the paper we use ACDM cosmo- 
logical parameters f2 m = 0.27, tt A = 0.73, Q b h 2 = 0.0225, 
h = 0.704, erg = 0.81, n s = 0.96, consistent with recent 
studies of the cosmic microwave background (Komatsu ct al. 
2011). 



2 SEMI-NUMERICAL REIONISATION 
MODEL 

Our model of inhomogeneous hydrogen reionisation simu- 
lates both the evolving ionisation state of the IGM and its 
thermal history. The former is calculated using the semi- 



Geil & Wyithe 


2008 


1 (see also 


Zahn et al. 2007 


Mesinger & 


Furlanetto|2007 


|Choudhury et al.|2009; Thomas et al. 2009 



for similar approaches). The ionisation state derived from 
this model is then used to evaluate both the ionising flux 
and the mean free path at the Lyman Limit, and to infer a 
temperature distribution for the IGM. 

2.1 Ionisation model 

We simulate the IGM inside a periodic, cubic grid of side- 
length L — 100/i _1 Mpc, containing N = 256 voxels per 
side-length, with an over-density field S(x) generated using 
the transfer function of Eisenstein & Hu (19991. The cal- 



culations detailing the evolution of ionisation structure in 
this density field (which is initialised to a neutral state at 
z = 99) are described in|Geil fc Wyithe| p008||. Briefly, the 



model evaluates the ionisation fraction Qs r ,r of a spherical 
region of scale R around a given voxel at redshift z as 



dQs r 



dt 



Aion 

0.76 



Qs R ,R 



O 



dt 

dF col (8 R ,R, 2,M min ) 



1 + S F 



dt 

D 1 (z) ' 



-Dl(Zobs) 



(l + z) 3 Qs B , R ,(l) 



where 8r is the region's over-density, ob is the case-B re- 
combination coefficient and C is the clumping factor. In this 
work we shall assume C = 3, broadly consistent with recent 



results from hydrodynamical simulations ( Pawlik et al. 2009 
Raicevic & Theuns 20111. The justification for this choice, 



as will be discussed in more detail in Section [5.1| is that it 
tends to provide the weakest overall constraints on reion- 
isation redshift. Meanwhile, Ai on , the number of ionising 
photons entering the IGM per baryon in galaxies, is left as 
a free parameter, the modelling of which will be described 
in more detail in Section [3] 

The production rate of ionising photons is then as- 
sumed to be proportional to Ni on and the rate of change 
of the collapsed fraction, Fed, in haloes above some min- 
imum threshold mass for star formation (M m i n in neutral 
and Mi on in ionised regions, corresponding t o virial temper- 
atures of 10 4 K and 10 5 K respectively, e.g. Dijkstra et al. 
|2004[ ). The collapsed fraction in a region of comoving ra- 
dius R and mean overdensity S(z) — 8{z h a )D\{z) / Di(z b B ) 
is found using the extended Press-Schechter model (Bond 
|et al.|199l] ), which gives 



F co i(5r, R, z) = erfc 



5c — 5 R (z) 



(2) 



numerical model developed by Wyithe & Loeb (20071 and 



where erfc(:r) is the complimentary error function, a\ is the 
variance of the overdensity field smoothed on a scale R and 
a"g al is the variance of the overdensity field smoothed on a 
scale -Rgai, corresponding to a mass scale of M m j n or Mi on 
(both evaluated at redshift z). 

Equation jl]) is solved at every voxel on a range of scales 
from L to L/N at logarithmic intervals of width AR/R = 
0.1. A voxel is then deemed to be ionised if it lies inside a 
region for which Qs r ,r > 1 on any scale R. By repeating this 
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procedure at successive redshifts, we obtain our ionisation 
field as a function of z. 



2.2 The ionising emissivity 

In order to compute the heating of the IGM, we must next 
obtain the ionising emissivity in each voxel of our simula- 
tion. We simulate the ionisation field at 30 regular redshift 
intervals between an initial redshift, Zo, and the reionisa- 
tion redshift, z r . These are, respectively, the redshifts corre- 
sponding to Q — 0.01 and Q — 1.0, where Q — Qs r ,r for 
the case where equation |TJ is solved at mean density and 
R — > oo. At each interval, if a voxel is flagged as ionised, we 
compute the intrinsic emissivity as the product of 7Vi on and 
the collapsed fraction above a minimum halo mass, Mj on . 
Using this, the proper ionising photon production rate (iV 7 ) 
in units of photons cm~ 3 s -1 , at each voxel, is 



iV 7 = 



p c n h (l + z) 3 (l + S x )N ion dF col (S x ,\m,z,M ion ) 



pm p 



dt 



(3) 



where p c is the present day critical density, m p is a proton 
mass and p is the mean molecular weight. In this instance, 
the smoothing scale A used to calculate the collapsed fraction 
is given by the minimum of the ionised bubble radius and 
the mean free path at the Lyman limit, Am- This is related 
to the ionising emissivity e„, (in units of erg s _1 cm -3 Hz -1 ) 
above the Lyman Limit frequency, vm, by 



AS 



dv, 



where we adopt a simple power law spectrum 



v 



em {vm < v < J4mi), 

(mell < V). 



(4) 



(5) 



Here eni is the emissivity at the hydrogen Lyman limit. Be- 
low the Hen ionisation edge, we take a spectral index of 
a = 3 characteristic of soft, population-II stellar sources 
(e.g. |Hui fc Haiman||2003"l |Furlanetto fc Oh|[2009l ). Above 
z/Heii, the spectrum is zero in order to enforce the require- 
ment that there be no significant stellar ionisation of He n. 
This is necessary for consistency with our semi-numerical 
reionisation model which only traces the propagation of hy- 
drogen ionisation fronts, and so cannot predict the larger, 
extended Hem fronts resulting from harder Hen ionising 
photons. On the other hand this is not an unrealistic re- 
quirement, since the double reionisation of helium is thought 
to not complete until z ~ 3 when the number density of 



quasars is at its peak (Madau & Meiksin 1994 Wyithe & 



Loeb 2003 Furlanetto & Oh 20081. Nevertheless, this as- 



sumption rules out the consideration of more exotic ion- 
isation models involving early populations of quasars or 
population-Ill sources, which could produce different ther- 
mal histories. More detailed modelling, which accounts for 
doubly ionised helium may therefore be performed in the 
future. In this work, however, we consider hydrogen reioni- 
sation by population-II stellar sources only. 



2.3 The ionising background specific intensity 

Having computed the intrinsic emissivity, the specific inten- 
sity, J„ , at each voxel may then be calculated by accounting 



for the propagation and attenuation of radiation. This is 
a complicated calculation since, in general, it depends on 
the topology of H n regions around each point in the simu- 
lation and the distribution of Lyman limit systems. A full 
numerical evaluation therefore requires radiative transfer, 
which is not included in our model. Fortunately, the long- 
term thermal evolution (for timescales longer than the local 
photo-ionisation timescale) is relatively insensitive to the ex- 
act amplitude of J„, and instead depends primarily on the 



spectral shape of the ionising radiation ( Hui & Haiman 2003 



Bolton et al. 2009). Therefore, to estimate the specific in- 



tensity we make the simplifying assumption that the mean 
free path is much smaller than the horizon scale, and use 



the local source approximation (e.g. Faucher-Giguere et al 



2009 1. We compute the mean free path in two regimes, prior 
to and following the overlap of Hn regions. 

Prior to overlap, the mean free path is set by the size 
of H n bubbles ( Gnedin fc Fan|2006 1 , and may therefore be 
found by partitioning the ionisation field into approximately 
spherical Hn regions. To do this, we smooth over the ion- 
isation field at successive scales, employing the same loga- 
rithmic descent as for our density field. At each scale we 
calculate the size of discontiguous regions with a standard 
friend-of-friends (FoF) algorithm. The bubble size for each 
voxel is then the scale of the smallest FoF region of which 
it is part, augmented by the smoothing scale. This has the 
effect of dividing bubbles that only barely overlap into sep- 
arate spheres so that the spherical approximation remains 
valid to a filling factor of around Q « 0.5. We then assume 
that the optical depth to ionising photons is infinite beyond 
the H n region radius. Lyman Limit photons therefore travel 
one mean free path, corresponding to the Hn region size, 
before being absorbed. 

In the post-overlap case, the IGM is fully ionised and 
the mean free path of Hi ionising photons is constrained 
by the abundance of Lyman limit systems rather than the 
sizes of Hn regions. In this instance, we employ a fit to the 
mean-free-path at the Lyman limit in proper Mpc, taken 



from Songaila & Cowie (20101 



Ai 



.6 



l + z 



r(2-j9) V 4.5 



Mpc, 



(6) 



where F denotes the gamma function, and /3 is the power law 
exponent of the Hi column density distribution. Observa- 
tionally, the value of this exponent is unconstrained at z > 6, 
however at lower redshift measurements indicate /3 w 1.5 
j Petitjean et al.|1993 Hu et al.|1995 1. Since this is roughly 
constant over the range 2 < z < 4 ( |Misawa et al.|[2007[ ), 
we simply extrapolate to high redshift and adopt /3 = 1.5 
throughout. We note that the column density distribution 
will change significantly as we enter the reionisation epoch, 
however as we are primarily interested in the value of /3 in 
the post-overlap regime, this does not affect our results. Fi- 
nally, we assume that a given H n region is in the pre-overlap 
regime if its radius is less than Am (i.e. the mean free path is 
not set by Lyman limit systems) and is in the post-overlap 
regime otherwise. 

Adopting the local source approximation, the ionising 
background specific intensity in each voxel may then be ap- 
proximated by 



4tt 



e v Ai/, 



© 0000 RAS, MNRAS 000, 000-000 



4 S. Raskutti et al. 



where \ v = Ahi(^/^hi) 3 '' 3_1 ' ) - Prior to overlap Am is fixed 
to the Hn region size, while post overlap it is given by 
equation Lastly, note that in both regimes at frequen- 
cies above the Lyman Limit, we have J„ oc l / _c "+ 3 (' 3 ~ 1 ). 
This effectively hardens the intrinsic stellar spectrum by 
a — > a — 1.5 , accounting for the additional heating due 
to filtering of the ionising radiation by discrete, Poisson dis- 



tributed absorbers in the intervening IGM (e.g. Miralda- 
Escude|2003) |Faucher-Giguere et al.|2009| >. 



2.4 The thermal evolution 

Once the specific intensity at a given voxel is known, the 
photo-ionisation rates, Li [s _1 ], and photo-heating rates, 
<7i [ergs -1 ], for i—H i, He i and He n, may be calculated 
via: 



T; = 



Si = 



(8) 
(9) 



where cn[y) and i>i are the photoionisation cross-section and 
frequency at the ionisation edge for species i, respectively. 
The temperature evolution is then computed by solving: 

( 7 -lW [G{t) _ A(t y 

kbP 



dT _ 
~dt ~ 



where y = |, and G — ^\ nigi and A are the total photo- 
heating and cooling rates per unit volume, respectively. We 
ignore any heating resulting from the growth of density fluc- 
tuations as this is a minimal effect throughout the short 



timescale of reionisation. Our code for solving equation ( 10 \ 
follows photo-ionisation and heating, collisional ionisation, 
radiative cooling, Compton cooling and adiabatic cooling 
for six species (H i, H n, He i, He n, He in, e~). The non- 
equilibrium abundances of these species are found by ad- 
ditionally solving three further differential equations cou- 



pled to equation ( 10 1 (e.g. Anninos et al. 1997 Bolton & 
Haehnelt 20071, and taking a uniform clumping factor of 
C = 3 for all species. We use the rates compiled in |Bolton] 
& Haehnelt ( 2007 1 with the exception of the case-B recom- 
bination and cooling rates of Hui & Gnedin ( 1997[) an d the 



photo-ionisation cross-sections of Verner et al. (19961. The 
gas density, p, is evaluated from our linear over-density field 
using a standard log-normal fit ( Coles fc Jones|1991[ ), which 
approximates the real non-linear density for the majority of 
voxels. Whilst this does not adequately reproduce the ex- 
treme over-dense and underdense tails of the density distri- 
bution ( Becker et al.|2007 1, its use does not represent a sig- 
nificant source of error since the temperature measurements 
we compare to are smoothed over the measured quasar near- 
zones (~ 5 proper Mpc). Inaccuracies in evaluating the ther- 
mal history of the small number of regions that fall within 
these tails are therefore unimportant. 

Using this prescription, we evolve the temperature of 
each voxel, independent of all others from an initial temper- 
ature To = 2.27(1 + z)K and an initially neutral state. This 
procedure yields a final temperature distribution over our 
simulation volume as a function of redshift. 



Table 1. Photo-ionisation rate constraints from the Lya forest 
opacity obtained by |Bolton et al. | ( |2005^ and |Wyithe fc Bolton] 
( |2011| l. The most likely derived value and errors at 68 per cent 
confidence intervals for the corresponding number of ionising pho- 
tons entering the IGM per baryon in galaxies, N- lon , are also given. 



Redshift 


r_i2 


1 Mon 


4 


n 7 _|_0.4S 


35.3±^2!i 


5 


0.47±g;3 


22.7±i 4 7 3 


6 


o.i8±°:o8 


13.0±« 5 ° 



2.5 Inclusion of quasar Hen photo- heating 

The final part of our model is including the effect of the 
quasars themselves on the thermal state of the IGM in their 
vicinity. The ionising flux from a quasar will significantly al- 
ter the thermal state of the surrounding IGM, as its harder, 
non-thermal spectrum will also reio nise He n (jBolton et al. 



20101. For any given quasar from Bolton et al. (20111, 
placed in the simulation volume at a random position, we 
account for this heating by augmenting the specific inten- 
sity in each voxel within the quasar proximity zone by an 
additional term 



exp 



R 
Am 



v N -3(8-1) 
VYU 



(11) 



where R is the proper distance from the quasar and L„. q = 
^Hi, q (^/^Hi) _aq is the quasar ionising luminosity, again 
taken to be a power law spectrum with index a q . There is 
considerable uncertainty concerning the exact value of this 
exponent towards high redshift. We therefore allow it to vary 
over the range 1.0 < a q < 2.0, which is consistent with 
recent indirect estimates from observations of quasar near- 



zone sizes at z > 6 ( Wyithe & Bolton 2011 1. To compute the 



quasar luminosities at the Lyman limit, we take the absolute 
AB magnitudes from Carilli et al. (|2010|). Finally, we fix the 



extent of the region which experiences additional ionisation 
and heating by each quasar to 5 proper Mpc, matching the 
radius within which the near-zones were analysed in [Bolton] 
et al. 1 12011 1. 



3 MODELLING OF IONISING SOURCES 

We next discuss the parametrisation of ionising sources in 
our model, encapsulated by the redshift-dependent parame- 
ter Mon, the number of ionising photons entering the IGM 
per baryon in galaxies. We begin by describing existing ob- 
servational constraints on the IGM ionisation state and then 
proceed to discuss how these are used to calibrate A r i on and 
its evolution with redshift. 



3.1 Observational constraints on A?i on 

3.1.1 The ionisation rate from the Lya forest opacity 

Observations of the Lya forest in quasar absorption spectra 
can be used to infer the ionising background due to luminous 
sources. The measured quantity is the mean transmission of 
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Lya flux, (F), along the line of sight, or equivalently the 
effective optical depth r e a = —log (F) . With knowledge of 
the physical properties of the IGM, these can be converted 
into the background photoionisation rate per hydrogen atom 
(Thi = r_i2 x 10~ 12 s _1 ). Using full hydrodynamical simula- 



tions of the IGM, Bolton et al. ( 2005 \ and Wyithe fc Bolton 



(20111 have measured T_i2 from the Lya forest opacity at 
,2 = 4 — 6. Their derived values are shown in Table [I] These 
measurements may then be compared to our own simulated 
ionisation rate, which is calculated using equation (JsJ . When 
combined with equations (|3|-([7|), and evaluated for the case 
of a uniformly ionised IGM at z < 6 this yields constraints 
on 7Vi on , which are also shown in Table [l] 



3.1.2 Thomson optical depth to CMB photons 

The electron scattering optical depth provides an integrated 
constraint on the IGM ionisation state throughout the epoch 
of reionisation. We use results from the latest WMAP obser- 
vations, tcmb = 0.088 ± 0.015 ( |Komatsu et~aT|2011[ ). This 
is related to the ionisation history by 



TCMB = 



Z CMB J, 

c^Q(z)(l + fne(z))nn(z)a T dz, (12) 
az 



where or is the Thomson scattering cross-section, 71h is the 
hydrogen number density and /n e is a correction factor due 
to the presence of ionised helium. In our model, we assume 
that the ionisation of Hei traces that of Hi, whilst Hen is 
instantaneously ionised at z = 3. 



3.2 Parametrisation of Ni on 

Using these observations, we may now constrain possible 
parametrisations of Ni on in our model. There is consider- 
able uncertainty concerning the exact shape of iVi on , par- 
ticularly for z > 6, where there are no direct observations. 
Indeed, several studies ( |Pritchard et al . 2009; Mitra et al.| 
20111 have found that Ni on is effectively unconstrained in 
the range 6 < z < 11. To explore the sensitivity of our 



simulated temperatures to changes in the ionising sources, 
we therefore choose to employ two different models for the 
redshift evolution of JVi on . 



as an a posteriori constraint on z r alongside the predicted 
temperatures in Section [5] 



3.2.2 Cubic Model 

We also explore arbitrary ionisation histories that best 
match the tcmb constraint, by employing a model that is 
redshift-dependent prior to z r . To ensure enough freedom to 
match all constraints, we choose to follow the approach of 
Pritchard et al. ( 2009 I and model Ni ^ as a cubic function, 



allowing for a large variety of different evolutions: 



Mon (Z) = (N p - N ) 



( z - go \ 3 _ ( Z - ZQ \ 

I Az j 6 \ Az ) 



+JVo,(13) 



where we have parametrised by the midpoint (zo,No), and 
the peak (zq + ^,iVp). F° r each reionisation redshift, we 
fit this curve to the measured electron scattering optical 
depth and the ionisation rate at z — 4 — 6. As in |Pritchard| 
|et al.| ( [2009| ), we impose a prior on the comoving cmissivity 
of N~, < 10 s1 s -1 Mpc -3 to disallow unphysically large ionis- 
ing backgrounds. We stress that this model is not physically 
motivated and clearly overfits the limited data we use. How- 
ever, as we do not seek to exactly find M on but only to view 
the effect of different ionisation histories on simulated IGM 
temperatures, this does not significantly impact upon our 
results. 



3.2.3 Optimal N ion models 

Finally, the best fit evolution of Ni on and the filling factor 
Q for the cubic and linear models are shown in Figure [I] 
for z r = 6.5 and z r = 9.0. For both cases, the behaviour 
of the emissivity matches quite well to results derived else- 



where (Choudhury & Ferrara 2006 Pritchard et al. 2009 



Mitra et al. 20111, although note that for early reionisation 



in the cubic model, star formation switches on later and the 
filling factor develops more quickly in order to match the 
electron scattering constraint. Importantly, the two mod- 
els do provide very different ionisation histories at each z r 
and therefore should allow us to probe the sensitivity of our 
temperature predictions to realistic variations in the ionising 
background. 



3.2.1 Linear model 

The simplest model, and one often used in previous sim- 
ulations of reionisation (e.g. Haiman & Loeb 1997 Geil 



|fc Wyithe||2008"l l, assumes a single, redshift independent 
parametrisation of A?i on . However, this cannot be reconciled 
with measurements of the photo-ionisation rate from the 
Lya forest at z < 6, so we employ a slight variant on the 
model. We take Ni on to be a constant up to the reionisation 
redshift, which is found by solving equation |l]) for a universe 
that is fully ionised by the desired redshift. This effectively 
fixes the reionisation history for any given z r . Then, to en- 
sure consistency with observations of F_i2, we allow Ni on (z) 
to evolve linearly once reionisation is complete in order to 
best fit the data points at z < 6. This is broadly consistent 
with predicted increases in the escape fraction towards high 
redshift ( Wyithe et al. 2010). As part of our constraint on 
z r we impose the resulting electron scattering optical depth 



4 RESULTS OF THE SEMI-NUMERICAL 
MODEL 

4.1 The IGM thermal state following hydrogen 
reionisation 

We now examine the variation of IGM temperature with 
reionisation redshift using the semi-numerical model de- 
scribed in the previous two sections. For any given z r , and for 
each of the two different parametrisations of Ni on , we com- 
pute an optimal Nion evolution that ensures Q — 1.0 by z r . 
The resulting ionising background then allows us to calculate 
the temperature field over the redshift range 5.5 < z < 6.5, 
when the quasars analysed in Bolton et al. (2011 1 are ex- 



pected to have had their most recent optically bright phase. 

The IGM temperatures predicted by our semi-numerical 
model at z = 6.0 are summarised in Figure [2] The top 
two rows display the results using our linear iVj on model, 



© 0000 RAS, MNRAS 000, 000-000 



6 S. Raskutti et al. 





z z 

Figure 1. Evolution of the number of ionising photons per baryon in galaxies which escape into the IGM, Ni an (left hand columns), 
and the Hn filling factor, Q (right hand columns) for reionisation histories with z r = 6.5 (top) and z r = 9.0 (bottom). Solid curves show 
the evolution for the linear A^on model, whilst the case of the cubic JVi on is shown by the dashed lines. The plotted points represent 
constraints imposed by values of T_i2 obtained from the Lya forest opacity by |Bolton et al~| ( |2005[ l and |Wyithe fc Bolton||2011[ |. The 
inset in the left hand panels gives the electron scattering optical depth for the linear and cubic N- lon models (r and t*, respectively). 



for reionisation completing by z r — 6.5 (first row) and 
z t — 9.0 (second row). The left hand panels display tem- 
perature maps through the centre of the lOOfo -1 Mpc simu- 
lation volume. The corresponding central panels show maps 
of the voxel ionisation redshifts relative to the redshift at 
which reionisation completes (i.e. Az — z — z r ). We note that 
the structure of ionised bubbles at fixed Az is very similar 
for the early and late models. In the case of late reionisa- 
tion, the temperature mimics this ionisation structure, with 
cooler regions centred on the small Hn bubbles that were 
ionised first, and have thus had more time to cool. However, 
regions of the IGM ionised during the overlap phase retain 
a clear imprint of heating. The situation is very different for 
early reionisation, where the temperature field shows little 
correlation with the reionisation topology. Instead, the tem- 
perature largely reflects the random Gaussian structure of 
the density field, indicating that even with reionisation oc- 
curring at z r — 9.0, the temperature has already settled to 
an almost asymptotic state by z = 6.0. 

This asymptotic limit is also clear from the 
temperature-density relation (right hand panels of FigurepJJ). 
For Zr = 9.0, the relation shows an almost power law depen- 



dence with 7 — 1 w 0.5 for T — To A 7 1 , characteristic of 
an asymptotic IGM thermal state which is set by the com- 



petition between photo- heating and adiabatic cooling (Hui 
|fc GnedTn"||1997| . However, in the case of late reionisation, 
the figure shows much wider 68 per cent confidence limits 
(dashed lines) below the highest densities, with an inverted 
power-law slope of 7 — 1 ~ —0.2 for low densities (A < 1.0). 
The scatter and anti-correlation between temperature and 
density, which arises due to the fact that lower density re- 
gions tend to be ionised last allowing them less time to cool, 
is in qualitative agreement with results from previous nu- 



merical and analytical studies ( Bolton et al. 2004 Trac et al. 
2008{ |Furlanetto fc Oh|2009| . 



The same set of results can also be seen for the cubic 
model of N lon in the two lower rows in Figure [2] In this 
case, the growth of ionised bubbles is significantly different 
at constant Az — z — z r for the two reionisation redshifts. 
As described in Section |3.2| late reionisation (z r = 6.5) is 
a more gradual process in this model, with ionisation times 
spread more evenly around the Hn bubbles. This has the 
effect of spreading out the IGM temperatures, so that re- 
gions ionised early tend to be cooler. Similarly, the scatter in 
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Figure 2. Maps of a two-dimensional slice through the centre of our 100/i -1 comoving Mpc simulation box at Z = 6.0 prior to the 
inclusion of any quasar heating. The upper six panels show results for the linear N- Ion model while the lower set show the same results 
allowing for a cubic Ni on . In each set, the upper row indicates maps for reionisation occurring at z T = 6.5 and the lower row for z T = 9.0. 
Left panels: The IGM temperature at the measurement redshift (z = 6.0). Central Panels: Maps of the ionisation redshift for each voxel 
relative to the redshift at which reionisation completes, Az = z — z r . Right panels: The IGM temperature, T, against overdensity, A, at 
z = 6.0, displaying both the median (solid) and 68 per cent intervals around the median (dashed). 
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the temperature-density relation is increased. On the other 
hand, for z r — 9.0, reionisation is far more rapid, effectively 
occurring over a redshift interval of Az — 3.0. Despite this 
difference, the temperature field evolves, by z — 6.0, to the 
same asymptotic state as for the linear TVion model (with 
a slightly reduced scatter in the temperature-density rela- 
tion due to the more abbreviated differences in ionisation 
time). This reaffirms the conclusion that the thermal mem- 
ory of hydrogen reionisation for z r — 9.0 is largely erased by 
z = 6.0. 

To see this more clearly, in Figure[3]we consider how the 
mean and standard deviation of the temperature distribu- 
tion at z — 6.0 evolve as a function of reionisation redshift. 
We note that the mean referred to here is not directly com- 



parable to the measurements of Bolton et al. (20111 as 



it is averaged over the entire density field in the simulation 
volume and is not simply the temperature at mean density. 
The results are shown for both the linear (solid red) and 
cubic (dashed blue) iVion models. The mean temperature 
predicted by the two models only varies by several hundred 
degrees Kelvin above z r w 9.0, and an asymptotic state is 
fully reached for z T « 10.5 — 11.0. However, for later reion- 
isation redshifts the temperature distributions vary quite 
rapidly, and for the linear A^ on model at z r = 6.0 the 1GM 
is ~ 5000 K hotter than the asymptotic value. This result is 
slightly more muted for the cubic model, as the elongated 
reionisation history means that large regions of the 1GM 
have longer to cool, thereby reducing the impact of heating 
on the mean temperature. However, even allowing for this, 
unit changes in the reionisation redshift correspond to tem- 
perature variations of up to 1500 K by z — 6.0. As this is 
comparable to the standard deviation of the distributions, 
changes of Az r w 1.0 should be distinguishable on the basis 
of volume averaged temperature measurements for z T < 9.0, 
but will require very precise data. 

Finally, note that the standard deviation exhibits a min- 
imum at z r ~ 7 for both models, but increases to both lower 
and higher redshift. At z r < 7, this scatter is driven by the 
large variation in temperatures resulting from recent, inho- 
mogeneous reionisation. In contrast, the scatter at z r > 7 is 
instead due to the steepening correlation between temper- 
ature and density which establishes itself following reion- 
isation. Both models asymptotically approach a value of 
a(T) ~ 3000 K for z r > 10. 

4.2 The temperature in quasar near-zones 



Bolton et al. (20111 used the Doppler parameters of 



Hi Lya absorption lines in the near-zones of 7 quasars at 
5.8 < z < 6.4 to measure the IGM temperature at mean 
density, To, within ~ 5 proper Mpc of these sources. In or- 
der to compare these observations to our model, we must 
therefore also include the effect of Hen photo-heating by 
the quasars themselves on the IGM temperature. 

We calculate this heating as follows. At any given 
voxel in our simulation volume, we switch on a quasar 
tq — 1.5 x 10 7 yrs prior to the redshift at which the quasar 
is observed, and use the procedure outlined in Section |2.5| 
to calculate the surrounding IGM temperature. We then 
draw 200 random lines-of-sight within a ~ 5 proper Mpc 
region around the quasar, and for each one we measure 
the temperature-density relation for all voxels. This relation 





Figure 3. Mean (top panel) and standard deviation (bottom 
panel) of the simulated temperature distributions at z = 6.0 prior 
to the inclusion of quasar heating. Note the results are obtained 
from the entire lOO/i" 1 comoving Mpc simulation box, and the 
mean is therefore volume averaged over all densities. We show 
results for two different cases; the linear JVi on model (solid red 
curves) and the cubic model allowing for an evolving Ni on (dashed 
blue curves). 



may then be used to infer a temperature at mean density 
averaged over all 200 lines-of-sight. Repeating this proce- 
dure for "quasars" centred on all voxels in the simulation 
box then gives a distribution for the temperatures around a 
single observed quasar. 

In Figure [4] we show the mean and standard deviation 
of the simulated temperature distributions at mean density 
for a quasar observed at z = 6.02. The results for the mean 
temperature are displayed for both hard (a q = 1.0, upper 
curves) and soft (a q = 2.0, lower curves) quasar spectral 
indices. The red solid and blue dashed curves again corre- 
spond to the linear and cubic Ni OI1 models. These retain all 
the features of the temperature distributions prior to quasar 
heating, including the asymptotic behaviour of the mean at 
z r > 9.0, and the minimum in the scatter at z T ~ 7.0. Note, 
however, the standard deviation tends to be slightly smaller, 
particularly at large z t . This is due to the fact that we are 
now only considering the average temperature at a single 
(mean) density. 
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Figure 4. Mean and standard deviation of the simulated distri- 
butions over temperature at mean density, p(Tq\z t ), in a quasar 
proximity zone at z = 6.02. Top: Simulated mean for three dif- 
ferent cases, the linear JVj on model (solid lines), the cubic model 
(dashed lines), and the cubic model accounting for the density 
bias of quasars (dot-dashed lines). We show simulated tempera- 
tures for o q = 2.0, which yield the lower curves and a q = 1.0 
giving the higher set of curves. Bottom: Standard deviations for 
the aforementioned models. In this case we only show results for 
oj q = 1.0, however standard deviations are similar regardless of 
quasar spectral index. 



We find that the maximum difference in the average 
temperature at mean density for the cubic iVi 0n model is 
only ~ 2 500K for early and late reionisation. Given that 
the 68 per cent measurement uncertainties on the tempera- 
ture around individual quasars obtained by |Bolton et al. | 
(2011[) are ~ 3 000-5 000K, it is unlikely that a strong con- 



straint can be obtained from any single line-of-sight observa- 
tion. However, the combination of all 7 quasars analysed by 



Bolton et al. (2011 \ reduces the 68 per cent confidence in- 



terval on the mean temperature of the quasar near-zones to 
~ 1 500K. Therefore, all seven lines-of-sight will be able, in 
principle, to distinguish between early and late reionisation 



scenarios. 



It is important to note that the rather uncertain quasar 
spectral index also introduces additional uncertainties into 
the interpretation of temperature. For example, in the cu- 
bic model the temperatures obtained for late reionisation 



combined with a soft quasar spectrum are almost identical 
to those obtained for early reionisation and a much harder 
quasar spectrum. Unfortunately, there are no direct obser- 
vations that constrain high-redshift quasar spectral indices, 
although a recent study of quasar near-zone sizes at z ~ 6 
demonstrated that the range 1.0 < a q < 2.0 is consistent 
with the data ( Wyithe fc Bolton|2011 ). This is in agreement 



with direct observations at lower redshift (e.g. Telfer et al 



|2002[ ), but in practice the quasar spectral index will remain 
an uncertainty in our analysis. 



4.3 Quasar Density Bias 

When calculating the measured temperature distribution ex- 
pected around the quasars, we have thus far assumed a ran- 
dom distribution of quasar positions in our simulation vol- 
ume. However, quasars and their host galaxies tend to trace 
high density regions owing to the clustering bias of massive 
haloes. Before calculating the constraints on z r , we therefore 
also model the quasar density bias. We weight the tempera- 
ture distribution calculated in Section [4. 2| by the likelihood 
jC s (5) of observing a quasar at that location 



p(T\w) 



]T 5(T(x)-T)£M*))- 



(14) 



The likelihood of observing a galaxy, and a corresponding 
resident quasar, is then estimated from the|Shcth & T ormen| 
(2002) mass function as 



{l + 5)u{l + v-'^) t 



-av*/2 



(15) 



where v = (<5 C - 8)/[a{M)], V = 8 c /[a(M)], and a = 0.707, 
p = 0.3 are constants. 

The results including this density bias within the cubic 
Mon model are shown by the dot-dashed curves in Figure 
g] There is a slight flattening of the temperature differences 
between different reionisation redshifts. This is due to the 
relatively early reionisation time of the high density regions 
where quasars are found, which allows them more time to 
cool to their asymptotic thermal state. At large values of 
z r the temperature is also systematically higher. This is be- 
cause the higher recombination rate in higher density regions 
results in slightly larger thermal energy input from photo- 
heating. This increased mean temperature makes it slightly 
more difficult for measurements of To to differentiate be- 
tween reionisation redshifts. However, this is a minimal ef- 
fect which will not strongly affect our constraints on z T . 



5 CONSTRAINTS ON THE REIONISATION 
HISTORY 

In the remainder of the paper we turn to calculating the con- 
straints that recent and possible future IGM temperature 
measurements at z ~ 6 impose on the hydrogen reionisation 
redshift. Aside from a number of fixed assumptions includ- 
ing the clumping factor C, star formation mass thresholds 
M m in,Mi on , the stellar spectral index a, and the evolution 
of iVi on , there are two critical free parameters in our model, 
2 r and a q . For any combination of these two free parame- 
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the i* near-zone temperature measurement, (^{z T , aV). Us- 
ing Bayes' Theorem, we then infer the likelihood for any set 
of 7 independent sightlines {q 1 }J^ 1 (corresponding to the 7 
measured temperatures from Bolton et al. |2011 |: 

p(T>\{q i }l =1 ,M)p({q i }U 1 \M) 



p({<Z l }J=i|D,M) = 



p(D|Af) 



(16) 



where D is the combination of observational constraints ob- 
tained for t cmb , the ionising background and the quasar 
near-zone temperature measurements, M is one of our mod- 
els for reionisation, and p(D\M) is a normalisation constant. 
The likelihood p(D|{g i }J =1 , M) can be easily found from 
the measured temperature distributions in Figure [3] and the 
measured errors on tcmb (which are taken to be Gaussian) 
and T_i2. The prior probability of choosing any set of sight 
lines p({q'}J =1 \M) is calculated from the density bias like- 
lihood in equation ( 15 1 and the assumption that prior in- 



formation on our free parameters is uniform over the ranges 
of interest, 1.0 < a q < 2.0 and 6.0 < z T < 12.0. Redshifts 
outside this range are not considered as they are effectively 



ruled out at the lower end by the Lya forest data ( Fan et al. 
2006[) and at th e upper end by the measurement of tcmb 
( Komatsu et al.||2011[). 

Using equation (16 I, we may then find absolute bounds 



on the reionisation redshift by summing over sets of possible 
sightlines 

p(z t \T>,M) = Y,P(W}i=i\'0,M)S[z I ({q i }J =1 \M)-z I ],(17) 

where the Kronecker delta 8 arises since all seven sightlines 
are drawn from a single simulation box and correspond to 
one value of z T . We note, however, that the same does not 
hold true for the quasar spectral indices, as each of the indi- 
vidual quasars may have different intrinsic spectra. There- 
fore, we only obtain constraints for each quasar spectral in- 
dex individually: 



pK|d,m) 



p({q i }I =1 |D,M)5K({q i }I =1 |M)-<](18) 



Equations ( 17 1 and ( 18 1 produce model dependent estimates 



for our parameters of interest. We note that there are sev- 
eral uncertain aspects of our modelling, which can lead to 
significant variations in sightline temperature, particularly 
for late reionisation (Figure Hb. 

To both compare the relative probability of our differ- 
ent models, and also obtain a less model dependent estimate 
of Zy- and a q , we use a Bayesian model comparison to calcu- 
late the posterior probability for model Mi given the data 
p(Ms|D) 



p(Mi|D) = 



p(D|M t )p(Afi) 
p(D) 



(19) 



Here p(D) is again a normalisation constant and we assume 
a priori that all models are equally likely so that p(Mi) is 
uniform. Model averaging can then be used to infer proba- 
bilities of free parameters across all models 



p( 2r |D) = ^p(z r |D,M i )p(M i |D), 



(20) 



which provides a more general set of constraints on z T and 
a q . These constraints are not completely model independent 



as the set of models we consider is not complete. However, as 
our models bracket most of the physical uncertainties, they 
should be reliable. 



5.1 Reionisation constraints from measured 
near-zone temperatures 

In the upper panels of Figure [5] we present the result- 
ing probability distribution functions over the parameter z T 
for the linear JV; on model. The distributions are computed 
for the two different sets of measurements presented by 
Bolton et al. ( |2011| , derived from their fiducial model (red 



solid curves) and a model with additional pressure (Jeans) 
smoothing in the IGM (blue dashed curves). The latter ef- 
fectively lowers the near-zone temperature constraints by 
~ 3000 K. The upper left panel displays the p(z T ) based on 
the measured values of To, tcmb and T_i2. For comparison, 
in the upper right panel we also plot the individual distribu- 
tions obtained by alternately considering constraints from 
either tcmb or To alone combined with our assumed Ni on 
evolution. 

As expected, the measurement of To alone does not 
differentiate between early reionisation redshifts, as for all 
scenarios the IGM has reached an asymptotic tempera- 
ture. However, the temperature data do provide constraining 
power for low z T . We find that irrespective of which set of 



temperature measurements from Bolton et al. ( 2011 1 is con- 



sidered, very late reionisation is ruled out to high confidence. 
This is because the 5 000K increase in IGM temperature for 
2 r < 6.5, pushes the mean quasar temperature well above 
the observations. 

In the lower panels of Figure [5] we show corresponding 
constraints for the cubic model excluding (lower left panel) 
and including (lower right panel) the effect of quasar den- 
sity bias. In this model the observed value of tcmb is used 
in defining the parameters governing A^ on . It is therefore a 
much weaker a-posteriori constraint and we can no longer 
rule out reionisation at high redshift. Similarly for low z T , the 
smaller temperature differences between reionisation scenar- 
ios means that late reionisation cannot be excluded with 
the same confidence. Furthermore, the relatively tight con- 
straints obtained when considering the temperatures derived 
assuming a larger Jeans scale are probably not robust, as the 
constraints from the fiducial measurements are favoured by 
a factor of around 3 over these. 

The final model averaged distribution is shown in Fig- 
ure [6] which corresponds to a mix of the fiducial results 
from the linear and cubic models. Even though this largely 
ignores the models with the tightest constraints, we may 
still conclude that z r > 7.9 (6.5) to 68 per cent (95 per 
cent) confidence. As a sanity check, we also verify that the 
optimal temperatures predicted by our reionisation models 
match the measurements quite well. Figure [7] shows the pre- 
dicted sightline temperature distributions for the maximum 
a-posteriori estimates of z r and a q compared to the observed 
quasar temperatures. We find that although two thirds of the 
modelled points lie within the 68 per cent confidence inter- 
vals of the observations, the model temperatures of quasars 
and 1 (following the numbering in Figure [7| appear to 
be systematically too high. However, for the fiducial |Bolton| 
|et al. | ( |2011[ ) measurements this does not provide a strong 
constraint ruling out our model, since the observed temper- 
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Figure 5. Constraints on reionisation redshift based on measured values of To, tcmb an d T-12 for different models of reionisation. 
In each case two probability distributions arc shown, representing constraints from the fiducial near-zone temperature measurements of 
|Bolton et al. | | |2011[ l (red solid curves) and allowing for a larger Jeans scale (blue dashed curves) . The relative probabilities of each 
model are shown in the top corner of each plot. Upper Left: Constraints from To and tcmb and the linear model for N Ion . Upper Right: 
Constraints obtained from the near-zone temperatures and N^ on alone. For comparison, the constraint from tcmb an d Mon alone is 
shown in black dot-dashed. Lower Left: The probability distributions for the cubic -/Vj on model. Lower Right: Distributions for the cubic 
model but now accounting for the bias of quasars to regions of high density. 



ature distributions are both positively skewed and broader 
than Gaussian. Therefore, although our model predictions 
lie on the edge of the 68 per cent confidence interval, they 
are only about a factor of 1.5 less probable than the peak. 
The same does not hold true for the second set of measure- 
ments which account for additional Jeans smoothing. These 
are ruled out with some significance, which may be indicate 
that the (conservative) systematic error due to additional 
Jeans smoothing adopted by |Bolton et al. | ([2011 ) may 
have been slightly overestimated. The only real outlier then 
comes from quasar 5 in the fiducial measurements, which is 
unsurprising since this provides a significantly higher mea- 
surement of proximity-zone temperature. This may indicate 
that a quasar spectral index which is harder than the prior 



we have assumed (i.e. a q < 1) may be required to reconcile 
this measurement with the simulations. 

It is interesting to contrast this result with recent evi- 
dence for a decline in the fraction of Lyman break galaxies 
which exhibit strong Lya emission from z = 5-7 ( |Ono et ah] 
2011||Pentericci et al.|2011||Schenker et al.|2011[ ). Assuming 
the Lycv emission escape fraction does not evolve rapidly, 
this decline may be attributed to a corresponding decrease 
in the ionised hydrogen fraction to Q « 0.5 by z = 7, imply- 
ing reionisation finished late. In addition, analyses of the size 
and shape of the near-zone observed in the spectrum of the 
quasar ULAS J 1120+0641 are consistent with Q < 0.9 at 
z ~ 7.1 ( |Mortlock et al.|201l[|Bolton et al.|2011[ ), assuming 
the red damping wing is not due to a proximate high column 
density absorber. Large neutral fractions of ~ 50 per cent 
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Figure 6. Final constraint on the rcionisation rcdshift based on 
measured values of To, tcmb and r_i2 when averaged over all 
our different models of reionisation. This includes different models 
for N lon , the bias of quasars to high density regions and the Jeans 
pressure smoothing scale. Model averaging significantly weakens 
constraints as it adds uncertainties to the temperature, but still 
rules out late reionisation to z r > 7.9 (6.5) at 68 per cent confi- 
dence and 95 per cent confidence respectively. 



at z ~ 7 may be difficult to reconcile with our temperature 
measurements, which indicate z T > 7.9 (6.5) to 68 per cent 
(95 per cent). Such large neutral fractions could perhaps be 
reconciled with our constraint if the average quasar spec- 
tral index is much softer than we have assumed, resulting in 
lower temperatures for a given z T due to less Hen photo- 
heating. On the other hand we note that other possible 
sources of hydrogen reionisation, such as population-Ill stars 
or mini-quasars, have harder spectra than we consider here. 
These would therefore tend to increase the heating imme- 
diately following reionisation. This would then drive IGM 



temperatures too high to match the Bolton et al. (20111 



observations assuming reionisation completed late. Conse- 
quently, sources which have harder spectra than we have 
assumed would strengthen rather than weaken our conclu- 
sion that late reionisation is unlikely. 

Similar considerations also apply to uncertainties in the 
clumping factor C, briefly mentioned in Section |2.1| This 
is a highly uncertain parameter, whose impact we did not 
fully model but instead simply assumed to take on a uniform 
value of C — 3. The clumping factor has an impact on tem- 
perature evolution, since increases in clumping factor tend to 
increase overall recombination, implying more photons are 
needed to achieve reionisation. This results in a correspond- 
ing increase in IGM heating on the order of ~ 1 000K for 
C = 5. However, we note that our choice of clumping factor 
lies at the low end of possible values derived from simula- 
Raicevic fc Theuns||201l"| ). Larger 



tion (Pawlik et al. 2009 



5.2 Prospects for constraints from future 
temperature measurements 

In the previous section we presented constraints on the 
reionisation redshift from the existing measurements of the 
IGM temperature around quasars at z ~ 6 (Bolton et al. 



20111. This work demonstrates the potential of high red- 



shift IGM temperature measurements to aid in ruling out 
very late reionisation, so it is therefore interesting to con- 
sider if future progress can be made with further data. 

We begin by creating mock catalogues of possible future 
observations. To do this, we assume that for a given z r , ionis- 
ing source model and density bias model, iV bs observations 
are to be drawn. For each one, the quasar spectral index 
is drawn from the uniform prior over a q and the quasars 
are assumed to have an optically bright phase which starts 
at a redshift drawn from uniform distribution in the range 
5.8 < z < 6.5. The median of each new observation is taken 
at random from our predicted temperature distribution for 
the chosen z T , the redshift at which the quasar turns on and 
Q q . Finally, the errors in the measurement around this me- 
dian are assumed to correspond to one of the 14 (7 fiducial 
and 7 augmented Jeans scale) observational measurements 
chosen at random. Once 7V b s measurements are drawn in 
this way for each z r and model, we then use equations ( 16 1 
- (1201) to obtain predicted probability distributions. 



values of the clumping factor should then only rule out late 
reionisation with greater significance for models which have 
the same z r . 



This procedure for the generation of mock observations 
is repeated 100 times, and the resulting distributions are 
averaged. We then marginalise over the different ionising 
source models to obtain a predicted distribution p(z T ) for ev- 
ery z r from which the measurements are drawn. This allows 
us to estimate the confidence with which the reionisation 
redshift can be measured as a function of the true z T . For 
simplicity, we omit the constraints imposed by tcmb in order 
to assess the new information drawn from the temperature 
observations alone. 

We present our results in Figure [8] The left panel shows 
the mode (solid line) and confidence intervals (filled regions) 
for -/V bs = 7, mimicking the constraints of the previous sec- 
tion. Meanwhile, the right panel shows the same results for a 
large number of observation^] iV b s = 100. We note that, in 
both cases, the mode typically recovers the true reionisation 
redshift in the models well, particularly for small z T . How- 
ever, at no reionisation redshift do we obtain a particularly 
tight constraint for seven sight-lines only. For late reionisa- 
tion, this is primarily due to uncertainties in modelling of 
both the quasar spectrum and the ionising efficiency. Simi- 
larly, if reionisation happened early, we will not be able to 
distinguish well between different z t as the IGM temperature 
reaches an asymptotic state independent of when reionisa- 
tion occurred. However, increasing the number of observa- 
tions from 7 to 100 can tighten our constraints on late reioni- 
sation considerably, increasing the potential 95 per cent limit 
on z r from 6.5 to 8.0. 

Improved constraints the average quasar EUV spectra 
index from low redshift data could help to tighten our prior 
on a q , and further improve limits at least on late reionisa- 



1 Note that a total of 100 quasar high-resolution spectra at z ~ 6 
is a very optimistic scenario. Obtaining one R ~ 40 000 spectrum, 
even for the brightest known quasars at these redshifts, requires 
~ 10 hours of integration on an 8-m class telescope. 
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Figure 7. Model averaged quasar temperature distributions at the maximum a-posteriori estimate of reionisation redshift. We show 
both the model distributions (blue) and the observed sightlines (red) with 68 per cent confidence error bars. In both cases, we also 
include the parameter estimates of z T and ct^ (rightmost point) through (leftmost point). For clarity, we omit redshift ranges and 
offset observations 3 and 5 by Az = 0.1. Left: Results obtained for the fiducial near-zone temperature measurements. Right: Results for 
the measurements with additional Jeans smoothing. 




Model z f Model z f 

Figure 8. Bias and uncertainty in the constraints on the reionisation redshift from near-zone temperatures, as a function of the true 
reionisation redshift in our models. In each case the mode of the predicted posterior distribution is shown (black line) alongside the 68 
per cent (lined region) and 95 per cent (hatched region) confidence intervals. Left: Results assuming 7 temperature measurements Right: 
Results for 100 near-zone temperature measurements, illustrating the optimistic case of large numbers of quasar sight-lines. 



tion. However, to tighten constraints on early reionisation 
using IGM temperature measurements ultimately requires 
observations which push to even higher redshifts. 



6 CONCLUSIONS 

Current understanding of the reionisation history is pri- 
marily built on existing observational evidence comprised 
of measurements of the Gunn-Peterson optical depth at 
z < 6, and the electron scattering optical depth of cosmic 



microwave photons. These measurements respectively probe 
only the tail end of the epoch of reionisation, and an inte- 
gral measure of the reionisation history. The available con- 
straints therefore permit a wide range of possible histories. 
Observations of the IGM temperature provide an additional 
complementary probe since the thermal memory of the IGM 
yields information at redshifts beyond the those where obser- 
vations are made. Measurements of the temperature around 
high redshift quasars therefore offer a renewed opportunity 
to constrain the reionisation history of the Universe. 

In this paper, we have used measurements of IGM tem- 
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perature at z ~ 6 (Bolton et al. 20111 to constrain the 



redshift at which the hydrogen reionisation epoch ended. 
Our analysis employs a semi-numerical model for reionisa- 
tion and heating, which is constrained to fit the CMB elec- 
tron scattering optical depth, and the ionisation rate inferred 
from the Lya forest at z — 6. We find that all models which 
reionise the Universe prior to z r m 9.0, settle to approxi- 
mately the same IGM thermal state by z = 6.0 with rela- 
tively little scatter in the temperature at mean density. For 
late reionisation, (z T < 8.0), the situation is different since 
there is significantly more variance in the IGM temperature 
distribution. This arises partly from differences in ionisation 
time, which allow for different cooling times and a greater 
scatter. However, we find that there is always a detectable 
thermal imprint of at least 2 500K, which differentiates the 
mean temperature signifying late reionisation z T « 6 from 
the asymptotic state achieved for z r > 9.0. This implies that 
measurements of temperature in quasar near-zones can be 
used to constrain the reionisation redshift if it is sufficiently 
late. 

Including observational constraints from the CMB elec- 
tron scattering optical depth and the Lya forest opacity, we 
conclude that z r > 7.9 (6.5) with 68 (95) per cent confi- 
dence. The z ~ 6 temperature measurements therefore rule 
out a very late completion to reionisation. The inclusion of 
further temperature measurements in the future will tighten 
these constraints. If reionisation occurred early, we estimate 
a sample of 100 quasar spectra would rule out z T < 8.0 to 
95 per cent confidence. However, in order to push these con- 
straints further, more detailed modelling of the IGM thermal 
evolution will also be required. In particular, the evolution 
of the ionising efficiency with redshift will have to be un- 
derstood in order to narrow the modelling uncertainties for 
later reionisation scenarios. Including the effect of hydro- 
gen reionisation sources which have harder spectra than we 
have assumed here will also be valuable. Finally, improved 
estimates of quasar spectral indices will aid in reducing un- 
certainties in the modelling of the near-zone temperatures. 
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